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ABSTRACT 

The motion of test particles in the gravitational fields generated by the f i rst fo ur 
members of the infinite family of generalized Kalnajs discs (| Gonzalez fc Reinal ((20 06)), 
is studied. In first instance, we analyze the stability of circular orbits under radial and 
vertical perturbations and describe the behavior of general equatorial orbits and so 
we find that radial stability and vertical instability dominate such disc models. Then 
we study bounded axially symmetric orbits by using the Poincare surfaces of section 
and Lyapunov characteristic numbers and find chaos in the case of disc-crossing orbits 
and completely regular motion in other cases. 
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1 INTRODUCTION 

A fact usually assumed in astrophysics is that the main part 
of th e mass of a typical spi r al gal axy is concentrated in a thin 
disk (|Binnev fc Tremainel (|l987l )). Accordingly, the study of 
the gravitational potential generated by an idealized thin 
disk is a problem of great astrophysical relevance and so, 
through the years, different approaches has been used t o ob- 
tain suc h kind of th in disk models(see Kuzmin (|l956l ) and 
Toomre (| 19631 ; 1 1964 ). as examples). A simple method to ob- 
tain the surface density, the gravitational potential and the 
rotation curve of thin disks of finite radius was developed by 
lHunterll| 19631 ). the simplest exa mple of d i sk obt ained by this 
method being the Kalnajs disk (jKalnais! (| 19721 )). which can 
also be obtained b y flat te ning a u n iform ly rotating spheroid 
i Wvse fc Mavalil j 19421) ; iBrandtJ (Il960l): iBrandt fc Beltonl 
( 19621 )). In a previous paper ([Gonzalez fc Reinal 12006)) we 
use the Hunter method in order to obtain an infinite family 
of thin disks of finite radius, an infinite family of generalized 
Kalnajs discs with a well behaved surface mass density. 

Closely related with the above study is the analysis of 
the motion of test particles in the gravitational field gener- 
ated by such disklike distributions of matter. In particular, 
the study of orbits in the equatorial plane is of clear astro- 
physical relevance due to its relation with the dynamics of 
intergalactic stellar motion or the flow of particles in accre- 
tion disks around black holes. Also, a knowledge of the disc's 
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internal kinematics is relevant for its subsequent statistical 
analysis, and the study of external particles motion help us 
to understand the behavior of stars belonging to the galaxy's 
remaining component (for example, the halo). Many of them 
cross back and forth through the disc, experiencing a fairly 
abrupt change in the gravitational force field. This fact gives 
rise to a large variety of chaoti c and regular orbits, as it was 
pointed out by Hunter l|2005l), i n the case of Kuz min-like 
discs, and Martinet et al (|l97ll ; 1 19731; 1 1974 Il975l ), in the 
case of Schmidt's models (Schmidt Il956l h 



In agreement with the above considerations, in this pa- 
per we shall focus on the kinematics around th e gene ralized 
Kalnajs discs, introduced by Gonzalez & Reina (|2006i ) . They 
form an infinite family of axially symmetric finite thin discs, 
whose first member is precisely the well-known Kalnajs disc 
i| 19721 ). The paper is organized as follows. First, in section 
[2] we present a summary of the main aspects of the gener- 
alized Kalnajs discs, the surface densities, the gravitational 
potential and the motion equations. Then, in section [3] we 
shall focus on equatorial orbits, i.e. z = trajectories (£ = 0, 
r\ = \J\ — R 2 /a 2 inside the disc and £ = \/R 2 /a 2 — 1, 77 = 
outside). At first we shall study the stability under radial 
and vertical perturbations of circular orbits, then examining 
the principal features and conditions of general equatorial 
orbits. In the next section, section [4] we present numerical 
solutions of ([?} describing some representative Poincare sur- 
faces of section. As it is expected, we find chaotic sections 
for disk-crossing orbits and regular sections for other cases. 
Some meridional plane orbits are plotted and the Lyapunov 
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characteristic numbers (LCN) are calculated. Finally, in sec- 
tion [5] we summarize our main results. 



2 THE GENERALIZED KALNAJS DISCS 

In this paper we shall focus on the kinematics around gen- 
eralized Kalnajs discs, introduced by Gonzalez & Reina 
(2006). They form an infinite family of axially symmetric 
finite thin discs, w hose fi rst member is precisely the well- 
known Kalnajs disc (I1972T ). The mass surface density of each 
model (labeled with the positive integer m) is given by 



£ m (7?) 



(2m + l)M 
2^2 



R 2 



2-i m-1/2 



(1) 



where M is the total mass and a the disc radius. Such 
mass distribution generates an axially symmetric gravita- 
tional potential, that can be written in terms of Legendre 
polynomials P n and second kind Legendre functions Q n as 
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We restrict our attention to these four members. The kine- 
matics corresponding to the remaining models (m > 5) is 
easily inferred from the features characterizing m = 1, 2, 3, 4. 
In all the calculations concerning with such models we 
choose a = M = G = 1, without loss of generality. 

Since each $ m is static and axially symmetric, the spe- 
cific energy E and the specific axial angular momentum £ are 
conserved along the particle motion. This fact restricts such 
motion to a three dimensional subspace of the (R, z, Vr, V z ) 
phase space. By defining an effective potential $ m as 



Here, — 1 < r\ < 1 and < £ < oo are spheroidal oblate co- 
ordinates, related to the usual cylindrical coordinates (R, z) 
through the relations 



R 2 = a 2 (l + £ 2 )(l-r 7 2 ) 



(3) 



The constants C*2n appearing in ([2]) are given by 
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Q2n (i£) an d G is the gravitational con- 



stant. Now, due to the presence of the term (m — n)! at the 
denominator, all the C2n constants vanish for n > m. 

According to ([2]), the gravitational potentials corre- 
sponding to the first four members are given by 
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the motio n wil l be determined by the equations (Binney & 
Tremaine [l987l ) 



R = V R 



z = V z , 



Vr = 
together with 



' dR ' 
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3 EQUATORIAL ORBITS 

The equilibrium points of the autonomous system (0 are 
Vr = V z = Z = 0, R = R c , where R c must satisfy the 
equation 



dR 



(Rc,0) 



' m \ or 



0. 



(9) 



(Rc,0) 



that is the condition for a circular orbit in the plane z = 0. 
In other words, the equilibrium points of ([7J occur when the 
test particle describes equatorial circular orbits of radius _R C , 
specific axial angular momentum given by 



£c = ±jR 3 c 



( 



V 9RJ (R 



(Bc,0) 



and specific energy 



E = $ m {R c ,0), 



(10) 



(11) 



where the subscript c in £ c indicates that we are dealing with 
circular orbits. 

In order to study the stability of these orbits under 
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Figure 1. Behavior of re 2 for m = 1,2,3,4. Values of R c such 
that this function is below the dashed line, corresponds to circular 
orbits that are unstable under radial perturbations. 




Figure 2. Behavior of u 2 for m = 1,2,3,4. Values of R c such 
that this function is below the dashed line, corresponds to circular 
orbits that are unstable under vertical perturbations. 



small radial and vertical (z-direction) perturbations, we an- 
alyze the nature of quasi-circular orbits. They are charac- 
terized by an epicycle frequency k and a vertical frequency 
v, given by (Binney & Tremaine Il987f ) 



d 2 $ m 
OR 2 



dz 2 



(12) 



(flc.o) 



This means that by introducing (|10p in the second deriva- 
tives of $ m we obtain k 2 and v 2 as functions of R c . Val- 
ues of R c such that k 2 > and (or) v 2 > corresponds 
to stable circular orbits under small radial and (or) verti- 
cal perturbations, respectively. Otherwise we find unstable 
circular orbits. The case m = 1 presents radial stability in 
the range < R c < 1 (k 2 — 3tt) but is radially unstable 
when 1 < R c < 1.198. For m — 2 and m = 3 we find 
radially unstable circular orbits with radius in the ranges 
2^2/3 < Rc < 1.075 and 2/VE < R c < 1, respectively. In 
contrast, circular orbits for m = 4 are always stable under 
small radial perturbations. We conjecture that models with 
m > 5 are also radially stable. Figure [1] shows the behavior 
of k 2 as a function of R c for m = 1, 2, 3, 4. Figure [2] show- 
ing the behavior of v 2 , illustrates the stability under vertical 
perturbations. We find the following ranges of vertical insta- 
bility: < R c < 1 for m = 1 (k 2 = -3tt/2); < R c < 0.943 
for m = 2; < R c < 0.688 for m = 3; < R c < 0.604 for 
m — 4. We see that the range of vertical instability decreases 
with ra. 

General equatorial orbits are determined by (O to- 



gether with the conditions V z = z = z = and E — 
Vjf/2 + <I> m (i?, 0). The motion is restricted by the inequality 
E > <t m (R, 0) and, in particular, we find bounded motion in 
a range Ri < R < R2 if it contains at least one critical value 
where $ m is minimum and $ m (i?i,0) < E < $ m (7?2,0). 
Figure [3] shows the effective potential for m = 2 and 
£ = 1.242, near the disc edge. At energies (a), (b), (c), 
and (d) we have bounded orbits, and (e), (f) corresponds 
to unbounded motion. In figure [3] we present the resulting 
phase portrait where the two regions of bounded and un- 
bounded motion are divided by a separatrix curve (dashed 
line). Similar phase portraits can be performed for m = 1 
and 777, = 3 if we set i a according to (|10[) . for 1 < R c < 1.198 
and 2/y/5 < R c < 1, respectively . For m = 4 the effective 
potential does not present local maximums and its phase 
portrait will not have any separatrix curve. 

At a given energy E, determined by equation (lll|l once 
we fix R c , the maximum possible value of the specific angular 
momentum is £ c - Therefore, it is convenient to parameterize 
£ by means of the ratio k = £/£ c (see next section). It is 
useful to calculate the range of values for £ c such that $ m 
has a minimum in < R < 1, i.e. ensuring that bounded 
motion is always possible inside the disk. In this way, we 
establish the limiting values of the integrals of motion for 
which particles will never escape from the source. Using the 
relation (|10[) . we found the following ranges for the specific 
axial angular momentum: (a) m = 1, < \£ c \ < \/3ty/2; (b) 
m = 2, < |4| < 2^10^/9; (c) m = 3, < |4| < ^/21tt/50; 
(d) 777 = 4, < \£ c \ < 15^7^/64. 



4 DISC-CROSSING ORBITS 

In this section, we present numerical solutions of motion 
equations (0 corresponding to bounded orbits outside the 
equatorial plane (except when they cross the plane z = 0). 
For certain values of E and £, they are confined to regions 
that contain the disc and will cros s back and forth through 
it. As it was showed by Hunter |2005h , this fact usually 
gives rise to many chaotic orbits due to the discontinu- 
ity in the z-component of the gravitational field, produc- 
ing a fairly abrupt change in their curvatures. There is an 
important exceptional case of this behavior: the Kuzmin's 
disc, characterized by an integrable potential of the form 
$ = -GM[R 2 + (a+\z\) 2 Y 1/2 , with a > 0. However, the so- 
called Kuzmin-like potentials, characterized by &(e) where 
e = [R 2 + (a + Iz]) 2 ] 1 / 2 , are non-integrable and present 
the behavior mentioned above. Generalized Kalnajs mod- 
els present a very similar structure and we can expect an 
analogous dynamics. Each potential & m (£,r]) can be cast in 
a Kuzmin-like form if we take into account that, according 
to ©, £ = (R+ + R^)/2a and 77 = (R+ - R-)/2ia, where 
R+ = [R 2 + (z + ia) 2 ] 1/2 and R- = [R 2 + (z -ia) 2 } 1/2 . More- 
over, they are characterized by a z-deri vative discontinuity 
in the disc, given by (Gonzalez & Reina 2006) 



d$> r , 
dz 



( d$>„ 



=0+ 



V dz 



= 2jrGS m (B). (13) 



Despite the above relation makes the KAM theorem inappli- 
cable, we also found a large variety of regular disc-crossing 
orbits. 

In Figure [5] we plot the level contours of for m — 
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Figure 3. Effective potential for m = 2 and I = 1.242. Hori- 
zontal lines corresponds to the following specific energy values: 
(a) -0.330975 (minimum of 4 2 ), (b)-0.330900, (c) -0.330850, (d) 
-0.33079 (maximum of i> 2 ), (e) -0.330700 and (f) -0.330600. Min- 
imum and maximum of $2 occur at R c = 0.925 and R c = 0.960. 




Figure 4. Phase portrait for m = 2, corresponding to the same 
values of I and E showed in Figure \3\ Here (a) and (d) are sta- 
ble and unstable circular orbits, respectively. Curves (b) and (c) 
corresponds to stable bounded motion, while (b) and (f) describe 
unbounded motion. The dashed line represents a separatrix be- 
tween bounded and unbounded motion regions. 



1,2,3,4, corresponding to E = -1.245 and I = 0.2, i.e. 
k = 0.276,0.266,0.263 and 0.262, respectively. For this val- 
ues, the motion of the particle is confined to a region con- 
taining the disc. The corresponding z = surfaces of sec- 
tion are shown in Figs. [6]|9] exhibiting a variety of regular 
and chaotic trajectories. The Fig. [6] that corresponds to 
the values determining the contour (a) in Fig. [S] shows a 
large KAM curve enclosing an island chain and three sets 
of three rings. There is also a stochastic region with two is- 
land chains near section's edge. The large KAM curve, the 
island chain and the central set of rings are produced by box 
orbits, while the lateral rings as well as the last two island 
chains are formed by loop orbits. The dotted curve, result- 
ing from a banana boxlet periodic orbit, divide the regular 
and stochastic region. 

Fig. [7] exhibits similar features as Fig. [5] This surface 
of section corresponds to the values defining the contour 
(b) shown in Fig. [5] In this case we see a defined central 
region of box orbits (four central rings) and an enclosing 
chaotic zone that contains a variety of resonant islands of 
loop orbits. In this case the regions of box and loop orbits 
are clearly separated, in contrast with Fig.[6]where they are 
alternated. The surface of section corresponding to m = 3 
(contour (c) of Fig. [5]) is showed in Fig. [8] exhibiting a regu- 




0.4 0.6 



Figure 5. Level contours of (a) $1, (b) $2, (c) <3?3 and (d) $4, 
when E = -1.245 and I = 0.2. 



Ill 


LCAf(±0.0001) 


1 


0.0108 


2 


0.0110 


3 


0.0118 


4 


0.0120 



Table 1. Estimation of the largest LCN for initial conditions 
z = 10~ 10 , R = 0.681,Vr = 0.819 of a disc-crossing orbit at the 
chaotic region of Figs. 16101 



lar region composed by a central zone of banana boxlets and 
two resonant island chains of loops. In the chaotic region we 
see island chains again and three denser zones near the sec- 
tion's edge, formed by a loop orbit. Some of those meridional 
plane orbits are plotted in Fig. 1101 Finally, Fig. [9] shows the 
surface of section for m = 4, k = 0.262 and E = —1.245. We 
find a very prominent chaotic region with only two island 
chains and a small regular region of box orbits. 

The stochastic regions in Figs. are due the over- 
lapping of ma ny resonances caused by presence of the disc 
(Hunter I2OO5V One can see this fact clearly in the three 
denser zones near the section's edge of Fig. [5] where the res- 
onant islands are almost overlapped at E — —1.245. When 
the energy increases to —1.215, for example, the overlap- 
ping is complete and the trajectory turns to be irregular 
(Fig. Ill[) . In Fig. [9] we also note the footprint traced by the 
overlapping of three prominent central islands. 

In order to quantify the degree of instability of the orbits 
we calculate their largest Lyapunov characteristic numbers 
(LCN), defined as 



LCN 



lim 

A D -> 
t — * 00 



ln(A/A D ) 



(14) 



where A and A are the deviations of two orbits at times 
and t respectiv ely. W e obtain LCN using the procedure of 
Benettin et al (1976J). Thus by fixing the motion integrals 
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Figure 6. Surface of section for m = 1, k = 0.276 and E = 
— 1.245. We have a small chaotic region with two resonant island 
chains and a central non-destroyed tori zone. 
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Figure 8. Surface of section for m = 3, k = 0.263 and E = 
— 1.245. We see a prominent chaotic zone with island chains en- 
closing a regular region of box and loop orbits. The latter corre- 
sponds to the two central resonant islands. 



2 



-2 




0.2 0.4 0.6 0.8 
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Figure 7. Surface of section for m = 2, k = 0.266, E = —1.245. 
The stochastic region is larger than in Fig. [5] and the non- 
destroyed tori zone is entirely formed by box orbits. 

asB= -1.245 and i = 0.2, choosing A D ~ 10" 9 and t = 
10 7 , we estimate N corresponding to a typical chaotic disc- 
crossing orbit for the cases m = 1, 2, 3, 4 (Table[l|. We found 
that the degree of instability increases modestly with m. 

An interesting phenomenon occurs when the effective 
potential has saddle points outside the disc. This happens for 
the m = 1 and m — 2 discs, in whose case there are critical 
points in the ranges 1 < R c < 1.198 and 1 < R c < 1.075, 
respectively. Such equilibrium points are outside the disc 
but near its edge, so for certain values of E the contour of 
$ m will contain only a fraction of the disc and a z = 
empty region. Then we will find bounded disc-crossing and 
non-disc-crossing orbits. For example, we choose k — 0.673 
(corresponding to R c — 1.116) and E = —0.335, obtaining 
the contour (a) of Fig. [T2j The resulting surface of section 



Figure 9. Surface of section for rrt = 4, k = 0.262 and E = 
— 1.245. The chaotic region is larger than in the above figures and 
contain only two resonant island chains (the dotted curves at the 
extremum island are due to a periodic figure-of-eight orbit). The 
regular central zone is made entirely by box orbits. 

(Fig. II 3p presents a large chaotic region to the left and a 
small totally regular region to the right. Both are divided 
by the saddle point, so for initial conditions near its left side 
the particle is "trapped" in the stochastic zone. In contrast, 
for initial conditions near the right side of the saddle point, 
the motion is confined to a region of non destroyed tori. In 
Fig. [T?] we show a similar situation for m — 2. This surface 
of section corresponds to k = 0.991 (R c = 1.043) and E = 
—0.389 (contour (b) of Fig. [12]). We find a large zone of disc- 
crossing orbits with a small chaotic component enclosing an 
island chain and regular rings of banana boxlets. To the 
right of the saddle point there is a small regular region. For 
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Figure 10. Orbits in the meridional plane for some initial condi- 
tions in Figure|8] (a) Loop orbit at R = 0.739172, V R = 0.362539 
(middle island chain inside stochastic region); (b) Banana boxlct 
orbit at R = 0.421542, V R = 0.362539 (the largest central ring); 
(c) Pretzel shape loop orbit at R = 0.405726, V R = 0.573413 
(second island chain inside regular region); (d) A loop orbit at 
R = 0.283155, V R = 1.57197 that forms the three denser zones 
near the section's edge. The outer boundary is the zero-velocity 
curve. 




0.20.40.60.8 1 0.20.40.60.8 1 

R R 



Figure 11. Orbits with the same I and initial conditions as in 
Fig. [10] (d). In the case (a), where E = -1.225, the orbit is still 
regular, whereas it turns chaotic when E = —1.215 in case (b). 
The outer boundary is the zero-velocity curve. 

the cases m = 3 and m = 4 it is not possible to obtain 
sections with the above features since their corresponding 
effective potentials have not saddle points (the product k 2 v 2 
is positive in R > 1). 



5 CONCLUDING REMARKS 

One remarkable fact, suggested by the analysis performed 
in Sect. 2, is that the stability of circular orbits under ra- 
dial and vertical perturbations increases with the parameter 
m (for m > 2). More larger m corresponds to more stable 
models, when we deal with circular orbits. Another impor- 
tant fact is that the range of axial angular momentum for 
which we can find bounded motion inside the disc, decreases 
with m. Thus, loosely speaking, more and more stable mod- 
els have less and less possibilities to maintain their particles 



Figure 12. Level contours of (a) *i, E = -0.335, £ = 1.287; 
(b) $2, E = —0.389, I = 1.196. Both cases have two connected 
regions, one of them containing the disc and the other is a disc- 
free zone. 



inside the disc. These considerations have special relevance 
in the search of equilibrium distribution functions (they are 
E, ^-dependent) characterizing such galaxy models. 

The numerical calculations showed in Sect. 3 confirm 
the analysis made by Hunter about disc-crossing orbits. 
There exist a chaotic motion induced by the presence of 
the disk and, despite the current versions of KAM theorem 
do not apply (due to the discontinuous force field), there 
is also a significant range of regular orbits. Moreover, since 
in our case we deal with finite disc models, a distinction 
between disc-crossing and non-disc-crossing orbits is some- 
times necessary. In the latter case we did not find chaotic 
motion, even in extreme situations where there are saddle 
points outside the field source, and the disc and disc-free 
region are connected (Figures [TBI and [14")) . 

Although one would be tented to think that such hyper- 
bolic exterior points can induce chaos at disc-free regions, 
what really happens is that the stochastic motion tends to- 
ward disc regions and a completely regular motion is devel- 
oped outside there. These considerations have special rele- 
vance in galaxy models with thin disc plus halo components. 
Particles belonging to the halo component will follow a mo- 
tion with the features mentioned above and, as it was showed 
by some authors, this fact determines decisively the internal 
structure of suc h stellar systems (see Ostriker, Spitzer and 
Chevalier (ll972T l). 
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Figure 13. Surface of section for m = 1, & = 1.287 and E = 
—0.335 (contour (a) in Fig. 11211 . There is a prominent chaotic 
zone of disc-crossing orbits to the left of R= 1.116 (saddle point) 
and a regular region to the right. 
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Figure 14. Surface of section for m = 2, I = 1.196 and E = 
—0.389 (contour (b) in Fig. 1121 . The disc-crossing region has a 
narrow stochastic zone and a prominent KAM-curves region. We 
have again a small regular zone to the left of saddle point at 
R = 1.043. 



REFERENCES 

Benettin, C, Galgani, L. and Giorgilli A., 1972, Phys. Rev. 
A, 14, 2338. 

Binncy, J. and Tremaine, S., 1987, Galactic Dynamics. 

Princeton University Press, Princeton, N. J. 
Brandt, J. C, 1960, Ap. J., 131, 211 
Brandt, J. C. and Belton, M. J. S., 1962, Ap. J., 136, 352 
Gonzalez, G. and Reina, J. 2006, MNRAS, 371 (4), 1873- 

1876. 

Hunter, C, 1963, MNRAS, 126, 299 

Hunter, C, 2005, Ann. New York Acad. Sciences, 1045 (1), 



